### Attempt at replicating BOHN'S Analysis of LAPOP 2007/2006
### A data file is needed, with the LAPOP survey. This can be obtained from lapop

### This script recodes variables, analyses the data, and produces several tables
### that are included inthe paper and in the webappendix
### it also produces a .RData object that is used in _replic20022006.R to plot the results
### of different surveys. But this .RData file is also available directly in this repository

### Needs to be ran sequentially

#### THESE ARE THE MAIN QUESTIONS IN LAPOP REGARDING VOTE RECALL
#VB2: Did you vote in 2006?
#BRAVB3: Who did you vote for in the first round?
#BRAVB4: Who did you vote for in the second round?
#VOL1: Did you vote in 2002?
#VOL2: Did you remember who you voted for in first round 2002? (only if yes in VOL1)
#VOL2A: Who did you vote fore in round 1 2002?
#VOL3: Did you vote in the second round of 2002?
#VOL3A: Who did you vote for in second round of 2002?


library(foreign)
library(car)
library(survey)
library(mvtnorm)
library(Zelig)
library(MatchIt)
rm(list=ls(all=TRUE))

setwd("path to tolder")
lapop2006 <- read.dta("brazil_lapop06.dta")  ### READ IN THE LAPOP 2007 DATA FILE
cat("Total observations in the data set:",nrow(lapop2006),"\n")  

### Missing data
tmp <- lapop2006[,c("VB2","BRAVB3","BRAVB4","VOL1","VOL2","VOL2A","VOL3","VOL3A")]
print(as.matrix(apply(is.na(tmp),2,sum)))

### Basic SEtup:
lula022.raw <- lapop2006$VOL3A=="Lula" #ignores all the missingenss, non respostes, etc...
lula062.raw <- lapop2006$BRAVB4=="Lula"
d <- data.frame(lula022.raw,lula062.raw) #create a new data frame to hold the new data
d$bf <- lapop2006$BF1A=="Sim" #Whether the person receives BF
d$bf.all <- lapop2006$BF1A=="Sim"|lapop2006$BF1B=="Sim"|lapop2006$BF1C=="Sim"|lapop2006$BF1D=="Sim" #Whether the person receives BF, BA, CA...
d$gender <- lapop2006$Q1=="Homem"
d$age.raw <-  as.factor(recode(lapop2006$Q2,"NA=NA;16:24=1;25:34=2;35:44=3;45:59=4;else=5"))
d$age <- lapop2006$Q2 #continuous variable
d$educ.raw <- recode(iconv(lapop2006$VS2,from="LATIN1"),"
	'Analfabeto / sem instru��o'=1; 
    c('Primeiro ano do ensino m�dio/ Primeiro ano do 2� grau',
      'Segundo ano do ensino m�dio/ segundo ano do 2� grau',
       'Terceiro ano do 2� grau/ Ensino m�dio incompleto',
       'Terceiro ano do ensino m�dio/Segundo grau completo')=3;
    c('Iniciou a faculdade/universidade mas n�o se formou',
       'Gradua��o/faculdade','Mestrado ','Doutorado','Lato-sensu')=4;
       else=2") #Bohn's coding
d$educ <- as.numeric(as.character(recode(iconv(lapop2006$VS2,from="LATIN1"),
	"'Analfabeto / sem instru��o'=0; 
     'Primeiro ano do ensino fundamental'=1;
     'Segundo ano do ensino fundamental'=2;
     c('Terceiro ano do ensino fundamental',
     'Quarto ano do ensino fundamental  incompleto')=3;
     'Quarto ano do ensino fundamental /Prim�rio completo'=4;
     'Quinta s�rie/ 1� ano gin�sio'=5;
     'Sexta s�rie/2� ano gin�sio'=6;
     c('S�tima s�rie /3� ano gin�sio','Oitava s�rie incompleta/4� ano do gin�sio incompleto')=7;
      'Oitava s�rie/4� ano gin�sio/Primeiro grau completo'=8;
     'Primeiro ano do ensino m�dio/ Primeiro ano do 2� grau'=9;
     'Segundo ano do ensino m�dio/ segundo ano do 2� grau'=10;
     'Terceiro ano do 2� grau/ Ensino m�dio incompleto'=11;
     'Terceiro ano do ensino m�dio/Segundo grau completo'=11;
      'Iniciou a faculdade/universidade mas n�o se formou'=12;
      'Gradua��o/faculdade'=15;
      'Lato-sensu'=16;
      'Mestrado '=17;
      c('Doutorado','P�s-doutorado')=19"))) #OUr coding
d$income.raw <- recode(iconv(lapop2006$VS13,from="LATIN1"),"c('Sem renda','At� R$ 260,00','De R$ 260,01 at� R$ 520,00')=1; 
                                       c('De R$ 520,01 at� R$ 780,00','De R$ 780,01 at� R$ 1.300,00')=2;
                                       c('De R$ 1.300,01 at� R$ 2.600,00')=3;
                                       c('De R$ 2.600,01 at� R$ 5.200,00')=4;
                                       c('De R$ 5.200,01 at� R$ 7.800,00','Mais de R$ 7.800,00')=5;
                                       c('N�o respondeu','N�o sabe')=NA") #Bohn's coding
d$income <- recode(iconv(lapop2006$VS13,from="LATIN1"),"
	   c('Sem renda','At� R$ 260,00','De R$ 260,01 at� R$ 520,00')='00-02 SM'; 
       c('De R$ 520,01 at� R$ 780,00','De R$ 780,01 at� R$ 1.300,00')='02-05 SM';
       c('De R$ 1.300,01 at� R$ 2.600,00')='05-10 SM';
       c('De R$ 2.600,01 at� R$ 5.200,00')='10+ SM';
       c('De R$ 5.200,01 at� R$ 7.800,00','Mais de R$ 7.800,00')='10+ SM';
                                       c('N�o respondeu','N�o sabe')=NA") #Our coding
d$region <- recode(lapop2006$ESTRATOPRI,"'Norte'='N';'Nordeste'='NE';'Centro-Oeste'='CO';'Sul'='S';'Sudeste'='_SE'")

### Now, performing recoding of vote variable as defined in the paper
d$lula062 <- recode(iconv(lapop2006$BRAVB4,from="LATIN1"),"'Nenhum (foi votar, deixou c�dula em branco/anulou voto)'='branco/nulo';
                                    c('N�o se lembra','NR')='missing';
                                    'Outro'='missing';
                                    'Inap (N�o votou)'='nao votou';
                                    c('Lula','PT')='lula';c('Geraldo Alckmin','PSDB')='alckmin'") #basic recoding
cat("Second round 2006:\n");print(100*prop.table(table(d$lula062)))
d$lula062 <- recode(d$lula062,"c('nao votou')=NA")
cat("Second round 2006 (excluding abstentions only):\n");print(100*prop.table(table(d$lula062)))
d$lula062 <- recode(d$lula062,"c('missing','nao votou','outro')=NA")
cat("Second round 2006 (excluding missing and abstentions):\n");print(100*prop.table(table(d$lula062)))
d$lula062.noblank <- ifelse(d$lula062=="lula",T,ifelse(d$lula062=="alckmin",F,NA))
d$lula062 <- d$lula062=="lula"


d$lula022 <- recode(iconv(lapop2006$VOL3A,from="LATIN1"),"NA='missing in file';
                                    'N�o votei'='nao votou';
                                    'Outro'='missing';
                                    'N�o lembro/ NS'='missing';
                                    'Branco/ Nulo'='branco/nulo';
                                    'Lula'='lula';'Serra'='serra'")
cat("Second round 2002:\n");print(100*prop.table(table(d$lula022)))
d$lula022 <- recode(d$lula022,"c('missing','nao votou','outro','missing in file')=NA")
cat("Second round 2002 (excluding missing and abstentions):\n");print(100*prop.table(table(d$lula022)))
d$lula022 <- d$lula022=="lula"


### Now, performing recoding of vote variables that are not used in regression, but reported in tables
d$lula061 <-  recode(iconv(lapop2006$BRAVB3,from="LATIN1"),"'Nenhum (foi votar, deixou c�dula em branco/anulou voto)'='branco/nulo';
                                              c('Lula  ( PT, PCdoB, PRB)','Lula','PT')='lula';
                                              c('Geraldo Alckmin (PSDB, PFL)','Geraldo Alckmin','PSDB')='alckmin';
                                              'Heloisa Helena (PSOL, PSTU, PCB)'='helena';
                                              'Cristovam Buarque (PDT)'='cristovam';
                                              c('NS/NR','N�o se lembra','NR')='missing';'Inap (N�o votou)'='nao votou';
                                              NA=NA;'Outro '='nsnr';else='outros'") 
d$lula061 <- factor(d$lula061,levels=c('branco/nulo','lula','alckmin','helena','cristovam','outros','missing','nao votou'))
d$lula061[iconv(lapop2006$VB2,from="LATIN1")=="Sim, votou"&is.na(d$lula061)] <- "missing"  ### These are people that say the voted, but did not declare in whom
d$lula061[iconv(lapop2006$VB2,from="LATIN1")=="N�o votou"] <- "nao votou"  ### These are people that say they didn't vote
d$lula061 <- recode(d$lula061,"c('missing','nao votou','missing in file')=NA")
d$lula061 <- d$lula061=="lula"

d$lula021 <-  recode(iconv(lapop2006$VOL2A,from="LATIN1"),"'Branco/ Nulo'='branco/nulo';
                                              c('Lula')='lula';
                                              c('Serra')='serra';
                                              'Ciro'='ciro';
                                              'Garotinho'='garotinho';
                                              c('NS/ NR','Outras respostas','N�o lembra','N�o lembro/ NS')='missing';c('N�o votei','Inap (N�o votou)')='nao votou';
                                              NA=NA;else='outros'")
d$lula021 <- factor(d$lula021,levels=c('branco/nulo','lula','serra','ciro','garotinho','outros','missing','nao votou'))
d$lula021[which(iconv(lapop2006$VOL1,from="LATIN1")
=="Sim, votou"& is.na(d$lula021))] <- "missing"  ### These are people that say the voted, but did not declare in whom
d$lula021[which(iconv(lapop2006$VOL1,from="LATIN1")=="N�o votou")] <- "nao votou"  ### These are people that say tthey did not vote
d$lula021[which(iconv(lapop2006$VOL1,from="LATIN1")=="Sim, votou"&lapop2006$VOL2=="Branco/ Nulo")]  <- "branco/nulo" #  #said they voted but recall having voted null or blank
d$lula021[which(iconv(lapop2006$VOL2,from="LATIN1")=="N�o votei")]  <- "nao votou" #  #said they voted but recall having voted null or blank
d$lula021 <- recode(d$lula021,"c('missing','nao votou','missing in file')=NA")
d$lula021 <- d$lula021=="lula"

cat("\n\nMissing Data in Set\n")
tmp <- d[,c("lula061","lula062","lula021","lula022")]
print(as.matrix(apply(is.na(tmp),2,sum)))


pre.data <- d[,c("bf","lula022","lula062","gender","age","age.raw","educ","educ.raw","income","income.raw","region")]
#pre.data <- subset(pre.data,is.na(lula2)==F)
print(as.matrix(apply(is.na(pre.data),2,sum)))



### REPLIC ATTEMPT #1: Analysis with raw (naive) dependent variable
tmp <- na.omit(d[,c("bf","lula022.raw","gender","age.raw","educ.raw","income.raw","region")]) 
cat("Raw set has",dim(tmp),"observations\n")
l1.raw <- zelig(lula022.raw ~   bf + gender + as.numeric(educ.raw) + as.numeric(income.raw) + as.numeric(age.raw) + region, model="probit",data=tmp)
tmp2 <- na.omit(d[,c("bf","lula062.raw","gender","age.raw","educ.raw","income.raw","region")]) 
cat("Raw set has",dim(tmp2),"observations\n")
l2.raw <- zelig(lula062.raw ~   bf + gender + as.numeric(educ.raw) + as.numeric(income.raw) + as.numeric(age.raw) + region, model="probit",data=tmp2)

### REPLIC ATTEMPT #2: Analysis with processed dependent variable but categorical independent variables
tmp <- na.omit(d[,c("bf","lula022","gender","age.raw","educ.raw","income.raw","region")]) 
cat("\n\n\n\nRegressions\n\n\n\n\n2002: Corrected DV set has",nrow(tmp),"observations\n")
l1.rep <- zelig(lula022 ~   bf + gender + as.numeric(age.raw) + as.numeric(income.raw) + as.numeric(educ.raw) +   region, model="probit",data=tmp)
print(summary(l1.rep))
l1.rep.alt <- zelig(lula022 ~   bf + gender + age.raw + income.raw + educ.raw +   region, model="probit",data=tmp)
print(summary(l1.rep.alt)) #enters all categorical variables as categorical variables
tmp2 <- na.omit(d[,c("bf.all","bf","lula062","gender","age.raw","educ.raw","income.raw","region")]) 
cat("2006: Corrected DV set has",nrow(tmp2),"observations\n")
l2.rep <- zelig(lula062 ~   bf + gender + as.numeric(age.raw) + as.numeric(income.raw) +as.numeric(educ.raw) + region, model="probit",data=tmp2)
print(summary(l2.rep))#enters all categorical variables as categorical variables
l3.rep <- zelig(lula062 ~   bf.all + gender + as.numeric(age.raw) + as.numeric(income.raw) +as.numeric(educ.raw) + region, model="probit",data=tmp2)
print(summary(l3.rep))
l2.rep.alt <- zelig(lula062 ~   bf + gender + age.raw + income.raw +educ.raw + region, model="probit",data=tmp2)
print(summary(l2.rep.alt))
tmp2.1st <- na.omit(d[,c("bf","lula061","gender","age.raw","educ.raw","income.raw","region")]) 
l2.rep.1st <- zelig(lula061 ~   bf + gender + as.numeric(age.raw) + as.numeric(income.raw) +as.numeric(educ.raw) + region, model="probit",data=tmp2.1st)
print(summary(l2.rep.alt))
tmp4 <- na.omit(d[,c("bf","lula062.noblank","gender","age.raw","educ.raw","income.raw","region")]) 
cat("2006: Corrected DV after discarding BLANK votes has",nrow(tmp4),"observations\n")
l4.rep <- zelig(lula062.noblank ~   bf + gender + as.numeric(age.raw) + as.numeric(income.raw) +as.numeric(educ.raw) + region, model="probit",data=tmp4)
print(summary(l4.rep))

### OUR RECODING OF CATEGORICAL VARAIBLES
tmp <- na.omit(d[,c("bf","lula022","gender","age","educ","income","region")]) 
cat("\n\n\n\nCorrected for Categorical Variables\n\n\n2002: Set has",nrow(tmp),"observations\n")
l1c <- zelig(lula022 ~   bf + gender + age + income + educ +   region, model="probit",data=tmp)
print(summary(l1c))
tmp2 <- na.omit(d[,c("bf","bf.all","lula062","gender","age","educ","income","region")]) 
cat("2006: Corrected DV set has",nrow(tmp2),"observations\n")
l2c <- zelig(lula062 ~   bf + gender + age + income + educ +   region, model="probit",data=tmp2)
print(summary(l2c))
l3c <- zelig(lula062 ~   bf.all + gender + age + income + educ +   region, model="probit",data=tmp2)
print(summary(l3c))
tmp2.1st <- na.omit(d[,c("bf","bf.all","lula061","gender","age","educ","income","region")]) 
l2c.1st <- zelig(lula061 ~   bf + gender + age + income + educ +   region, model="probit",data=tmp2.1st)
print(summary(l2c.1st))
l3c.1st <- zelig(lula061 ~   bf.all + gender + age + income + educ +   region, model="probit",data=tmp2.1st)
print(summary(l3c.1st))
tmp4 <- na.omit(d[,c("bf","lula062.noblank","gender","age","educ","income","region")]) 
cat("2006: Corrected DV with BFALL set has",nrow(tmp3),"observations\n")
l4c <- zelig(lula062.noblank ~   bf + gender + age + income + educ +   region, model="probit",data=tmp4)
print(summary(l4c))

######### MAtching
cat("\n\n\n\nMATCHING\n\n\n2002\n")
tmp <- na.omit(d[,c("bf","lula022","gender","age","educ","income","region")]) 
matching <- matchit(bf~region+income+gender+age+educ
                   ,data=tmp,exact=c("income","region"),replace=T,ratio=2)   
md <- match.data(matching)
md$income <- factor(md$income) #refactor because high bracket droped
l1m <- zelig(lula022 ~   bf + gender + age + income + educ +   region, model="probit",data=md)
print(summary(l1m))

cat("\n\n2006\n")
tmp2 <- na.omit(d[,c("bf","bf.all","lula062","gender","age","educ","income","region")]) 
matching <- matchit(bf~region+income+gender+age+educ
                   ,data=tmp2,exact=c("income","region"),replace=T,ratio=2)   
md <- match.data(matching)
md$income <- factor(md$income) #refactor because high bracket droped
l2m <- zelig(lula062 ~   bf + gender + age + income + educ +   region, model="probit",data=md)
print(summary(l2m))
matching <- matchit(bf.all~region+income+gender+age+educ
                   ,data=tmp2,exact=c("income","region"),replace=T,ratio=2)   
md <- match.data(matching)
md$income <- factor(md$income)
l3m <- zelig(lula062 ~   bf.all + gender + age + income + educ +   region, model="probit",data=md)
print(summary(l3m))

bf0 <- setx(l2m,bf=F,region="_SE",gender=F,income="00-02 SM")
bf1 <- setx(l2m,bf=T,region="_SE",gender=F,income="00-02 SM")
fd.m <- summary(sim(l2m,x=bf0,x1=bf1))$qi.stats$rr

tmp2.1st <- na.omit(d[,c("bf","lula061","gender","age","educ","income","region")]) 
matching.1st <- matchit(bf~region+income+gender+age+educ
                   ,data=tmp2.1st,exact=c("income","region"),replace=T,ratio=2)   
md.1st <- match.data(matching.1st)
md.1st$income <- factor(md.1st$income) #refactor because high bracket droped
l2m.1st <- zelig(lula061 ~   bf + gender + age + income + educ +   region, model="probit",data=md.1st)

cat("\n\n2006 Matching droping blank votes \n")
tmp4 <- na.omit(d[,c("bf","lula062.noblank","gender","age","educ","income","region")]) 
matching <- matchit(bf~region+income+gender+age+educ
                   ,data=tmp4,exact=c("income","region"),replace=T,ratio=2)   
md <- match.data(matching)
md$income <- factor(md$income)
l4m <- zelig(lula062.noblank ~   bf + gender + age + income + educ +   region, model="probit",data=md)
print(summary(l4m))


###### ASSEMBLE THE TABLE FOR PAPER
stack.zelig <- function(x,name="reg"){                                         #x is a summary(zelig()) object                         
 core<-stack(data.frame(round(rbind(x$coefficients[,"Estimate"],x$coefficients[,"Std. Error"],x$coefficients[,"Pr(>|z|)"]),3)))
 core$ind<-as.character(core$ind)
 names(core)[1]<-name
 core$type <- rep(1:3,nrow(core)/3)
 core<-rbind(core,c(length(x$deviance.resid),"zN",NA))
 return(core)}
 
nice.table <- function(x){
  x <- x[c(grep("bf",x$ind),grep("gender",x$ind),grep("age",x$ind),
            grep("income",x$ind),
            grep("educ",x$ind),grep("region",x$ind),grep("Intercept",x$ind),
            grep("zN",x$ind)),]
    x$type <- recode(x$type,"1='';2='se';3='p-val';NA='N'")
    x$ind <- ifelse(x$type=="",x$ind,x$type)
    x <- subset(x,ind!="se")
    x <- x[,-2]  
    return(x)
}
 
t1<-merge(stack.zelig(summary(l1.rep),"r2002"),stack.zelig(summary(l2.rep),"r2006"),by=c("ind","type"),all=T)
t2<-merge(stack.zelig(summary(l1c),"c2002"),stack.zelig(summary(l2c),"c2006"),by=c("ind","type"),all=T)
t3<-merge(stack.zelig(summary(l1m),"m2002"),stack.zelig(summary(l2m),"m2006"),by=c("ind","type"),all=T)
the.table <- merge(t1,merge(t2,t3,by=c("ind","type"),all=T),by=c("ind","type"),all=T)
the.table <- nice.table(the.table)
write.csv(the.table,file="table-regressions.csv", na = "",row.names=F)

##### ASSEMBLE ADDITIONAL RESULTS FOR 2006 FOR APPENDIX
##### With Bohn's specificaition and categorical variables and dropping blanks
ta1 <- merge(stack.zelig(summary(l1.rep.alt),"r2002alt"),stack.zelig(summary(l2.rep.alt),"r2006alt")) #alternative
ta2 <- merge(merge( #dropping blanks
        stack.zelig(summary(l4.rep),"r2006noblank"),stack.zelig(summary(l4c),"c2006noblank"),by=c("ind","type"),all=T),
        stack.zelig(summary(l4m),"m2006noblank"),by=c("ind","type"),all=T)

ta3 <- merge(merge( #first round
        stack.zelig(summary(l2.rep.1st),"r20061st"),stack.zelig(summary(l2c.1st),"c20061st"),by=c("ind","type"),all=T),
        stack.zelig(summary(l2m.1st),"m20061st"),by=c("ind","type"),all=T)
            
ta4 <-  merge(merge( #expanded BF.ALL
        stack.zelig(summary(l3.rep),"r20061bfall"),stack.zelig(summary(l3c),"c20061bfall"),by=c("ind","type"),all=T),
        stack.zelig(summary(l3m),"r20061bfall"),by=c("ind","type"),all=T)
        
table.a1 <- nice.table(ta1)
table.a2 <- nice.table(ta2)
table.a3 <- nice.table(ta3)
table.a4 <- nice.table(ta4)
write.csv(table.a1,file="table-appendixregressions01.csv", na = "",row.names=F)
write.csv(table.a2,file="table-appendixregressions02.csv", na = "",row.names=F)
write.csv(table.a3,file="table-appendixregressions03.csv", na = "",row.names=F)
write.csv(table.a4,file="table-appendixregressions04.csv", na = "",row.names=F)


#### VOte in 2002 and Vote in 2006 for comparison with other surveys
#### Stores predicted values in a format that can be added to figures #### Make sure this 
vote.evs <- list()
d$income <- recode(iconv(lapop2006$VS13,from="LATIN1"),"
	   c('Sem renda','At� R$ 260,00')='00-01 SM';
	   c('De R$ 260,01 at� R$ 520,00')='01-02 SM'; 
       c('De R$ 520,01 at� R$ 780,00')='02-03 SM';
       c('De R$ 780,01 at� R$ 1.300,00')='03-05 SM';
       c('De R$ 1.300,01 at� R$ 2.600,00')='05-10 SM';
       c('De R$ 2.600,01 at� R$ 5.200,00')='10+ SM';
       c('De R$ 5.200,01 at� R$ 7.800,00','Mais de R$ 7.800,00')='10+ SM';
                                       c('N�o respondeu','N�o sabe')=NA")                                   
tmp <- d[,c("lula021","gender","age","educ","income","region")] 
tmp$income <- factor(tmp$income)
comp2002 <- zelig(lula021 ~   gender + age + educ  + region + income , model="logit",data=tmp)
evs.vote <- evs.vote.lb <- evs.vote.ub <- NULL
    for(k in 1:length(levels(tmp$income))){#do simulation and store values for categories of income
        x.out <- setx(comp2002, income=levels(tmp$income)[k], region="_SE", gender=F)
        evs <- try(summary(sim(comp2002, x = x.out))$qi.stats$ev,silent =TRUE) #caso haja colapso, nao deve ocorrer com menos categorias em id
        evs.vote <- rbind(evs.vote,evs[,1])
        evs.vote.lb <- rbind(evs.vote.lb,evs[,3])
        evs.vote.ub <- rbind(evs.vote.ub,evs[,4])
    }
    rownames(evs.vote) <- levels(tmp$income)
    evs <- list(vote=evs.vote,ub=evs.vote.ub,lb=evs.vote.lb)
    vote.evs[["LAPOP2007-2002"]]<-evs  


tmp <- d[,c("lula061","gender","age","educ","income","region")] 
tmp$income <- factor(tmp$income)
comp2006 <- zelig(lula061 ~   gender + age + educ  + region + income , model="logit",data=tmp)
evs.vote <- evs.vote.lb <- evs.vote.ub <- NULL
      for(k in 1:length(levels(tmp$income))){#do simulation and store values for categories of income
        x.out <- setx(comp2006, income=levels(tmp$income)[k], region="_SE", gender=F)
        evs <- try(summary(sim(comp2006, x = x.out))$qi.stats$ev,silent =TRUE) #caso haja colapso, nao deve ocorrer com menos categorias em id
        evs.vote <- rbind(evs.vote,evs[,1])
        evs.vote.lb <- rbind(evs.vote.lb,evs[,3])
        evs.vote.ub <- rbind(evs.vote.ub,evs[,4])
    }    
    rownames(evs.vote) <- levels(tmp$income)
    evs <- list(vote=evs.vote,ub=evs.vote.ub,lb=evs.vote.lb)
    vote.evs[["LAPOP2007-2006"]]<-evs    

save(vote.evs,file="predictedvoteLULA_LAPOP2007.Rdta")    









